Analytical framework for recurrence-network analysis of time series 
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Recurrence networks are a powerful nonhnear tool for time series analysis of complex dynamical 
systems. While there are already many successful applications ranging from medicine to paleoclima- 
tology, a solid theoretical foundation of the method has still been missing so far. Here, we interpret 
an £-recurrence network as a discrete subnetwork of a "continuous" graph with uncountably many 
vertices and edges corresponding to the system's attractor. This step allows us to show that various 
statistical measures commonly used in complex network analysis can be seen as discrete estimators 
of newly defined continuous measures of certain complex geometric properties of the attractor on 
the scale given by e. In particular, we introduce local measures such as the e-clustering coefficient, 
mesoscopic measures such as £-motif density, path-based measures such as e-betweennesses, and 
global measures such as e-efficiency. This new analytical basis for the so far heuristically motivated 
network measures also provides an objective criterion for the choice of e via a percolation threshold, 
and it shows that estimation can be improved by so-called node splitting invariant versions of the 
measures. We finally illustrate the framework for a number of archetypical chaotic attractors such as 
those of the Bernoulli and logistic maps, periodic and two-dimensional quasi-periodic motions, and 
for hyperballs and hypercubes, by deriving analytical expressions for the novel measures and com- 
paring them with data from numerical experiments. More generally, the theoretical framework put 
forward in this work describes random geometric graphs and other networks with spatial constraints 
which appear frequently in disciplines ranging from biology to climate science. 

PACS numbers: 05.45.Tp, 89.75.Hc, 05.45.-a 



I. INTRODUCTION 

Analogies are a fundamental motor of innovation in 
physics and other disciplines, since they allow the trans- 
fer of theoretical insights, results and techniques from 
one field to the other. In the last years, complex network 
theory has been particularly successful in providing uni- 
fying concepts and methods for understanding the struc- 
ture and dynamics of complex systems in many areas of 
science, ranging from power grids over social networks to 
neuronal networks [iHl]- Similarly, nonlinear time series 
analysis aims to gain insights on a wide variety of nat- 
ural, technological, and experimental dynamical systems 
drawing on a generic body of theory and methods ,5J . 

By exploiting analogies in the structure and descrip- 
tion of complex networks and dynamical systems, a num- 
ber of new network-based techniques for nonlinear time 
series analysis have been proposed recently [5]. The first 
class of these methods makes use of graph representations 
of certain similarity relationships between state vectors 
or groups of state vectors (e.g., cycles) in phase space. It 
includes transition networks based on a coarse-graining 
of phase space [7 , cycle networks [S], correlation net- 
works [9,, fc-nearest-neighbor [10], and adaptive nearest 
neighbor networks [TTl [T^ as well as e-recurrence net- 
works [131 [13] . The latter three techniques harness the 
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fundamental analogy between the Poincare recurrence 
structure [TS] of a time series in phase space, which is 
commonly represented by a binary recurrence matrix and 
allows us to recover basic dynamical invariants of the 
underlying system j.l6| , and the binary adjacency ma- 
trix describing a complex network. Other methods such 
as visibility graphs 17 work in the time domain and 
focus on studying stochastic properties of time series. 
Aside from these network-based approaches for investi- 
gating (possibly multivariate) time series from isolated 
dynamical systems, efforts have been spent for develop- 
ing techniques for studying fields of time series such as 
functional brain networks [T51 - BD] in the neurosciences or 
climate networks [^Tl - Bl] in climatology. In summary, all 
methods mentioned above propose a mapping from the 
time series to the network domain and then proceed to 
interpret the statistical properties of the resulting (usu- 
ally complex) network in terms of the underlying system's 
dynamical properties. 

While these interpretations are mostly based on empir- 
ical findings for paradigmatic model systems and heuris- 
tic arguments, only a few rigorous results are available. 
So far, Lacasa et al. have pointed out a relationship 
between the scaling exponent of the degree distribution 
Pk{k) ex k~^ in visibility graphs constructed from frac- 
tional Brownian motion and the Hurst exponent [25j . 
Furthermore, close relationships between the transitivity 
properties (network transitivity and local clustering co- 
efficients) [26] as well as the degree distribution's power- 
law scaling exponent 7 [27] of e-recurrence networks and 



the (fractal) global and local dimensionality of the at- 
tracting set underlying the time series have been found. 
Constituting random geometric graphs [55] , e-recurrence 
networks represent the geometry induced by the time se- 
ries in phase space in a simple and well-defined way. This 
enabled Donner et al. [Mj [2^ to define continuous tran- 
sitivity properties depending solely on the geometry of 
the (attracting) set 5* and the probability density func- 
tion p{x). These can in turn be calculated analytically for 
paradigmatic model systems with smooth and self-similar 
geometry and are approximated by the corresponding 
discrete e-recurrence network measures. Notably, most 
kinds of time series networks proposed so far are spatial 
networks [29] , since vertices are embedded either in phase 
space or on the time axis, implying that general results 
obtained for this class of networks are applicable to time 
series networks as well. 

Recently, e-recurrence networks have been demon- 
strated to be a particularly useful tool in diverse applica- 
tions of nonlinear time series analysis ranging from model 
systems [Q [El [H [Ml [30H36] via experimental data [STj- 
HD] to recent and paleo-climate records [3 [131 [HHSl as 
well as financial time series |31j . They allow us to un- 
cover complex bifurcation scenarios [131 E] and to re- 
liably distinguish between chaotic and non-chaotic dy- 
namics [33 . Furthermore, the local and global transi- 
tivity characteristics of e-recurrence networks have been 
shown to enable us to trace unstable periodic orbits [2] 
and to define alternative notions of fractal dimension |26j 
independently of earlier approaches. An important ad- 
vantage of nonlinear e-recurrence-network-based time se- 
ries analysis is that it performs well with significantly 
shorter time series (0(10^) data points [13l [33l [41] ) than 
required by classical techniques like estimating the max- 
imum Lyapunov exponent from data [SJ [31] • This ren- 
ders e-recurrence networks readily applicable to the anal- 
ysis of non-stationary real world data. The method has 
also been applied successfully to time series with irreg- 
ular sampling and/or uncertain timing of observations 
that are commonly found in the geosciences or in astro- 
physics [41, 42_. 

An e-recurrence network is completely defined by its 
adjacency matrix Aij{e) which is obtained from a (mul- 
tidimensional) time series x{ti), i — \, . . . ^N ,hy 



A,j{e) = Q{e-\\x{U)-x{tj)\\)-5^, 



(1) 



where Q{) is the Heaviside function, e a threshold used 
for defining the neighborhood of a state vector x{ti), || • | 
some norm, and 5ij denotes Kronecker's delta introduced 
to avoid self-loops in the network. Given univariate ob- 
servational or experimental time series, it is usually nec- 
essary to reconstruct the corresponding system's trajec- 
tory in some higher dimensional phase space to recover 
its recurrence structure reliably (e.g., by time-delay em- 
bedding OSB]). 

Within the recurrence network, vertices represent ob- 
servations or state vectors in phase space, while edges 
indicate a close proximity between two state vectors. 



Recurrence networks and their statistical properties are 
related, but complementary to the established concepts 
of recurrence plots (the recurrence matrix is given by 
Rij{s) = Aij(s) + Sij) and recurrence quantification anal- 
ysis (RQA) [16j. In contrast to RQA, which considers 
temporal dependencies between observations in form of 
diagonal and vertical line structures in the recurrence 
plot, recurrence network analysis discards all temporal 
information and solely quantifies the geometry of the un- 
derlying set S (e.g., an attractor) [l4ll4T] . 



Given the diverse and successful applications of e- 
recurrence network analysis reported in the literature, 
it is important to establish a firm theoretical founda- 
tion for advancing the understanding of the method. 
Building on earlier work [26], we propose here an an- 
alytical framework for e-recurrence network analysis of 
time series encompassing neighborhood-based transitiv- 
ity measures, mesoscopic measures relying on network 
motifs [47^ , path-based network characteristics as well as 
spectral and random walk-based measures. Specifically, 
our theory describes all graph-theoretical recurrence net- 
work quantifiers that have been used in the literature so 
far [M: . Beyond forming a solid theoretical basis for this 
modern nonlinear approach to time series analysis and 
fostering its detailed understanding in a way compara- 
ble to that of standard linear time series analysis [48], 
our analytical framework opens several avenues for prac- 
tically improving the method when dealing with finite 
(real-world) time series: (i) We are able to obtain closed- 
form analytical results for paradigmatic model systems 
with stochastic (uniform and Gaussian noise) and deter- 
ministic (periodic, quasi-periodic, and chaotic) dynamics. 
These can in turn be harnessed as a benchmark for the 
discrete standard estimators from complex network the- 
ory which have been employed so far [2], e.g., for assess- 
ing the estimators' bias and variance, (ii) This bottom- 
up approach allows us to design improved, weighted sta- 
tistical estimators }49j which may be more appropriate 
in specific situations, (iii) Moreover, our framework en- 
ables us to derive rigorous bounds for feasible values of 
the recurrence threshold e, the most important param- 
eter of the method, the choice of which is critical when 
analyzing finite (experimental) time series [3D]. We will 
argue in Sec.|V]that the concepts and measures developed 
in this paper can be readily generalized to describe the 
structure of a wider class of spatial networks, e.g., ran- 
dom networks with an arbitrary prescribed edge length 
distribution pi{l). 

This paper is organized as follows: We introduce a 
continuous framework for recurrence network analysis in 
Sec. jnl After reviewing the corresponding discrete esti- 
mators (Sec. Ill), we present examples ranging from peri- 
odic and quasi-periodic dynamics and higher dimensional 
symmetric sets over chaotic maps to stochastic processes 
and compare some of the results to discrete estimates 
(Sec. IV). We conclude with a discussion of these achieve- 
ments (Sec.fv]). 




FIG. 1: Illustration of a set 5* (gray), where g(x,y) denotes 
the geodesic distance between x,y £ S. 



are closely related to a certain notion of the fractal di- 
mension of an underlying set S and its associated prob- 
ability density p(x) with x G S |2S1- To capture this 
theoretically, continuous versions of both measures de- 
noted C{x;e) and T(e) have been defined together with 
a continuous degree density p{x; e). 



1. Local measures 



Definition 1. The continuous e-degree density 



II. CONTINUOUS FRAMEWORK 

A. General setting 

Let us consider a path-connected Lebesgue-measurable 
subset 5" C X of an rn-dimensional compact smooth man- 
ifold X with a non- vanishing continuous probability den- 
sity function p : S ^ (0,oo) with Jgd"^xp{x) = 1. 
We will use the abbreviation J dfi{x) = Jgd"^xp{x) 
throughout the paper, where /x is a probability measure 
on S (Fig. T]). Then we can define "continuous" equiv- 
alents of all relevant graph-theoretical measures for e- 
recurrence networks which may be approximated by cal- 
culating their discrete counterparts in the limit e ^ 0, 
TV — )■ oo (Sec. III). Here e is the threshold used for net- 
work construction (Eq. (fl])) and N denotes the number 
of data points (samples, phase space vectors, . . . ) consid- 
ered. These measures capture the properties of a "con- 
tinuous" network with uncountably many vertices and 
edges which may be defined by a continuous analog of 
the adjacency matrix, the adjacency function 



A{x, y) = e(e - ||a; - y||) - 5xy 



(2) 



for all x,y € S. It is important to realize that the 
framework introduced in this paper is not restricted to 
e-recurrence networks alone, but may be more generally 
applied to describe random geometric graphs (aka spatial 
networks) |251 [SDl [HI] and other types of networks with 
spatial constraints ^l [SU [S3] . 

In the following we will formally define the pro- 
posed continuous recurrence network measures and dis- 
cuss their properties, interrelationships and interpreta- 
tions (see Table ^ for an overview). Statements made for 
the limits £ — > and x -^ y for x,y (z S should be under- 
stood to hold for smooth S and p. We do not consider 
tliem for fractal geometries explicitly. 



B. Neighborhood-based measures 

Among other interesting properties, it has been shown 
recently that the local and global transitivity properties 
of e-recurrence networks measured by the local clustering 
coefficient Ci and the global transitivity T, respectively, 



p{x;e) 



dn{y) 



(3) 



B,(x) 



measures the probability that a point y randomly drawn 
according to p lies in an e -neighborhood B^{x) = {y G 
S ■■ ||a:-y|| < e} ofx. 

Definition 2. In turn, the continuous local e-clustering 
coefhcient of any point x G S, 

C x; e = '-^ -— ^ , 4 

is the probability that two points y and z randomly drawn 
according to p are closer than e given they are both closer 
than e to X. 

Definition 3. The continuous e-matching index 
Ai{x,y,e) measures the overlap of the neighborhoods of 
x,y G S, 



M{x,y;e) 



JbJx 



{x)nB,{y) 



dp(z) 



L 



B,(x)UB,{y) 



dp{z) 



(5) 



It gives the probability that a point z drawn randomly 
from Bi;{x) according to p is also contained in B^{y) and 
vice versa. 

For a; — > y, M{x, y; e) — > 1. Furthermore, M{x, y; s) = 
Oif||x-2/||>2e. 



2. Global measures 



Definition 4. The continuous e-edge density 



Pie) 



dp{x)p{x;e) 



= I dp{x) / dp,{y) 

Is JB^ix) 



(6) 



is the expectation value of the continuous e-degree density 
p{x;e). 



TABLE I: A summary of the continuous geometric properties of the set S and its associated probabihty density p defined in 
Sec. In] See main text for formal probabilistic interpretations. 



Class 

Neighborhood- based 
Local 



Global 



Mesoscoptc 



Path-based 
Local 



Global 



Name 

Continuous e-degree density (Eq. |3|) 
Continuous local e-clustering (Eq. |4|) 
Continuous e-matching index (Eg. pi) 
Continuous e-edge density (Eq.jGIl) 
Continuous e-transitivity (Eq. |7|) 
Continuous global e-clustering (Eq. Jsk) 
Continuous e-assortativity (Eq. 



Continuous e-motif density (Eq. (lOl) 



Continuous e-closeness (Eq. |l2|) 
Continuous e-efficiency (Eq. |l3|) ^^ 

Continuous e-shortest path betweenness (Eq. |l4r 
Continuous e-shortest path edge betweenness (Eq. 
Continuous e-average path length ( Eq. (^25|) 
Continuous glob al e- efficiency (Eq. (|27|) 
e-diameter (EqJ32 



Heuristic interpretation 

Local density 
Local dimension [26| 
Local density gradient 
Average local density 
Global dimension |26j 
Average local dimension 
Average local density gradient 



Higher-order density structure/ 
Density anisotropy 



r Geometric centrality 
rri, r Geometric bottleneckishness 



} 



Average separation 



e-radius (Eq. (33l) 



Geometric diameter 
Geometric radius 



Definition 5. As a global measure of geometric transi- 
tivity, we define the continuous e-transitivity of S as 

T{e) = dfi{x)dti{y)dfi{z)Q{e-\\x-y\\)x 



xe(£-||y-z||)e(£-||z-a;||) 



d^ix) dti{y) d^i{z)e{e - \\x - y\\) x 



xe(e-||z-a;||) , (7) 

which is the probability that among three points x,y,z 
drawn randomly according to p, y and z are closer than 
e given they are both closer than e to x. 

Definition 6. Similarly, the continuous global e- 
clustering coefficient 



C{e) 



d^{x)C{x;e) 



is the expectation value of the continuous local e- 
clustering coefficient C{x;e) (Eq. 01) j. 

Note that the above defined measures of transitivity 
have been mainly considered for the supremum norm Lq^ 
in [25]. 

Definition 7. Continuous e-assortativity 

A{e)=r(p{x;e),piy;e) \ \\x ~ y\\ < s] (9) 



gives the Pearson product-moment correlation coeffi- 
cient 148^ of the degree densities p{x; e) and p{y; e) of 
all points x,y that are closer than e to each other. 

A{£) can be considered as a measure of the smoothness 
of the set S and the probability density p [H] . In the limit 
e — ?► we have A{e) — >■ 1. 



C. Mesoscopic measures 

Motifs of order a are small connected subgraphs of 
a vertices that are embedded within the topology of a 
complex network |47| . For combinatorial reasons, usu- 
ally a < 5 is considered. Measuring motif densities is 
a useful approach for quantifying higher-order neighbor- 
hood relationships in complex geometries |llj and may 
be seen as a generalization of the transitivity concepts 
introduced above. 



\°) Definition 8. The continuous e-motif density 



M^ 



flfdfiixA J] e(e 



Xk\ 



(10) 



quantifies the frequency of occurrence of a certain recur- 
rence motif of order a described by the corresponding edge 
set E'a, where (3 = 1, . . . , n{a) and n{a) is the total num- 
ber of distinct motifs of order a. M^ is the probability 
that a points drawn randomly according to p are arranged 
according to the recurrence motif described by E^ . 



For example, the density Afy of the recurrence motif U 
of order 4 is measured by 



M^ 



II I dfi{xA e{e - \\xi - X2\\) X 

Xe(£ - \\X2 - X3\\)e{e - \\X3 - Xi\\). 



In contrast to the study of motifs for other general com- 
plex networks 47 , it is neither meaningful nor necessary 
to normalize motif densities by their expectation values 
for randomized networks here. The reason is that the 
Ma already have a natural probability interpretation. To 
render the results for different a more comparable, one 
may consider to use relative motif densities normalized 

We conjecture that motif densities as generalizations of 
the continuous e-transitivity are like the latter related to 
certain notions of the dimensionality of the set 5* and its 
associated probability density p [55] . This would allow us 
to define and study a new class of motif-based measures 
of dimensionality analogously to the sequence of Renyi 
dimensions from dynamical systems theory |54) . 



D. Path-based measures 



called a global geodesic or shortest curve [SS]. Depending 
on the geometry of 5", there may be none, one, or multi- 
ple distinct global geodesies connecting x and y, but in 
a sufficiently well-behaved set S, there will usually be a 
unique global geodesic for almost every pair x,y (Fig. fTl) 
and almost every metric (in particular, for the Euclidean 
metric), where by "almost every" we mean as usual that 
the set of exceptions has zero measure. Note that, how- 
ever, for some pathological metrics global geodesies are 
rarely unique, including the Li and Loo metrics. 

To understand the reasoning behind the following def- 
initions, one has to note that discrete shortest paths of 
the e-recurrence network approximate global geodesies 
connecting two points x,y G S for small e and large A'^. 
Then the shortest path length lij (e) (the minimum num- 
ber of edges that have to be traversed to reach vertex i 
from vertex j [1]) approximates g{x{ti),x{tj)), i.e.. 



£lij{e) sa g{x{ti),x{tj)) 



(11) 



(Fig. II]), where x 
£ -^ 0, and if N 



x{ti) and y = x{tj). In the limit 
> oo sufRciently fast, we argue in 



Appendix |A| that indeed slij{e) — >■ g{x{ti),x{tj)), inde- 
pendently of which metric is used for constructing the 
£-recurrence network. 



While the neighborhood-based properties defined 
above describe the small-scale geometry of the set S 
and probability density p, path-based measures quantify 
their global geometry in terms of global geodesies (see 
Fig. [2] for examples). Most of the path-based concepts 
defined below do not conceptually depend on the thresh- 
old e. Nevertheless, we introduce the appropriate scal- 
ing with £ into the definitions for consistency with the 
corresponding discrete estimators from complex network 
theory (Sec. Ill) and the recurrence network literature. 

As in standard topological terminology, a path in S is 
a continuous function / : [0, 1] — > 5, and its path length 
l{f) e [0,00] is the supremum of J2'i=id{f{U-i), f{ti)) 
over all n > and all tuples = to < ■ ■ ■ < tn = 1, where 
d(-,-) is some metric. Note that l{f) can be infinite in 
which case the path is called non-rectifiahle. For points 
X, y, the geodesic distance g{x, y) G [0, 00] is the infimum 
of l{f) over all paths in 5* from x to y (i.e., with /(O) — x 
and /(I) = y) ^ ■ A corresponding path of this length is 



^ Note that when defined in this way, g{x, y) may change discon- 
tinuously under continuous changes of the probabihty density p. 
This is because we require the geodesies to stay within S which 
consists of all points x where p(x) 7^ 0. When p{x) is continu- 
ously changed to zero, the length of geodesies running through x 
for p{x) > may change abruptly once p{x) = is reached, e.g., 
when X constitutes some kind of geometric bottleneck. If this be- 
havior is undesirable, one may consider generalized p-weighted 
notions of the geodesic distance. These could be motivated by 
an analogy to the optical path length in heterogeneous and non- 
isotropic media in physics, where the probability density p would 
play the role of the spatially varying refractive index. 



1. Local measures 

Definition 9. Given that a point y is drawn randomly 
according to p, continuous £-closeness 



c(x;£) 



9{x,y) 



dn{y) 

s £ 



(12) 



e I / dii{y)g{x,y) 

is the inverse expected geodesic distance of y to another 
chosen point x in units of e. 

Definition 10. Similarly, continuous local £-efRciency 
e{x;e) = |^dM2/)(^^) ' (13) 



= £ / dfi{y)g{x,y) 
Js 

gives the expected inverse geodesic distance of y to x mea- 
sured in units of e. 

Both c{x; e) and e{x] e) quantify the geometric close- 
ness of X to any other point in S given a probability den- 
sity p. Hence, points in the center of 5* will carry larger 
values of c{x] e) and e(x; £) than those on its boundaries 
(see Fig. [2] and below) . 

Definition 11. Continuous £-shortest path betweenness 

b{x- e) = ff dniy) dn{z) ^i^^^l^ , (14) 

JJs cr{y,z) 



c(x;e) 



e(x;s) 






bjx;s) 



U1 




FIG. 2; Local path-based measures continuous e-closeness c{x;e), -efficiency e{x;e), and -shortest path betweenness bi{x) 
(based on Eq. (15l), 62(2;; e) (based on Eq. (I61) in three example sets 5* with a uniform density p, two convex sets (circle and 
square), and a nonconvex set. Grayscale indicates the value of the measures (white: small, black: large) obtained by Monte 
Carlo numerical integration using the Euclidean norm for small e to avoid boundary effects. Note the more complex structure 
of the betweenness field, displaying particularly large values (dark) at the inward corners where many shortest paths must 
pass. In contrast to the path-based measures shown here, continuous e-degree p{x;e) and local e-clustering coefficient C{x;e) 
are constant in the interior of S due to the uniform p. Variations of these measures due to boundary effects occur only closer 
than £ to the boundary of S [26] . 



is the probability that a point x lies on a randomly chosen 
global geodesic connecting two points y, z drawn randomly 
from S according to p. Here, a(y, z;x;e) denotes the 
number of times x £ 5 lies on a global geodesic between 
y, z G S and cr(jj, z) is the total number of global geodesies 
between y, z (Fig. fj). 

In pathological situations, e.g., for certain open sets 
5", ^{y, z) may be zero even when the geodesic distance 
g{y, z) is well-defined and finite. We ignore these cases 
for now. 

There are several ways to formally define criy, z]X]e). 
Using a parametrization /^(i) of the family of global 
geodesies connecting y and z, with t e [0, 1] and /^(O) = 
y, /k(1) = z, we may write 



ai{y,z;x;e) = ai{y,z;x) 



^■(j/.z) 

E 

K = l 



dtSif^it)~x), (15) 



where S{-) is Dirac's multi-dimensional delta function. 
Alternatively, we can include the finite e-effect by 
counting all shortest paths that pass through the e- 
neighborhood of x by setting 



a2{y,z;x;e) 



E 



dtQ{e-\\f^{t)-x\\). (16) 



Both variants of (T(y, z; x] e) yield different, yet qualita- 
tively similar results for b{x; e) as is illustrated in Fig. (l2|. 

Given convex domains S, u{'y, z) = 1 always holds, 
i.e., there is only one straight line connecting y and z, 
parametrized by f[t) — y+t{z — y). For one-dimensional, 
convex sets S and using (Ti(y, z; x), continuous e-shortest 
path betweenness simplifies to 



b{x) ^2 Ij dn{y) dii{z)Q[x - y)Q{z ~ x). (17) 



Definition 12. Continuous e-shortest path edge be- 
tweenness 

h{x,y-e)^ //d/i(z)dMH ''^^'7''^'f^\ (18) 

JJs (j{z,w) 

is the probability that two points x, y both lie on a ran- 
domly chosen global geodesic connecting two points z, w 
drawn randomly according to p. a{z,w;x,y;e) counts 
the number of global geodesies between z, w which con- 
tain X, y. 



Analogously to continuous e-shortest path betweenness 



B 





FIG. 3: Illustration of the definition of continuous £-shortest 
path betweenness (the set S is indicated by gray shading). 
(A) There are a{y,z) = 2 global geodesies connecting y,z £ 
S, but only a{y, z; x;e) = 1 includes x £ S. (B) In this 
example, x lies on all four global geodesies between x and y, 
i.e., a{y,z) = a{y,z;x;£) = 4. 



It should be noted that for general 5* and p, the center 
of mass 



X 



dfi{x)x 



(23) 



does not necessarily extremize c(a;;e), e{x;s) or b{x). 
However, for convex S the generalized continuous e- 
closeness 



Cjj yX^ 6 ) 



dfi{y) 



9{x-,y) 



»7\-i 



(24) 



dy^{y)9{x.yT 



can be shown to assume a global maximum at x = X for 
the special case 77 = 2. In turn, the standard continuous 
£-closeness c(x;e) (77 = 1) is maximized at the geometric 
median or Fermat- Weber point j 57| . 



&(x;e), we can define this quantity as 
ai{z,w;x,y;e) = ai{z,w;x,y) 

<t{z,w) . I 

= E / dt5{f.{t)-x)] X 



X I j^ dt6iUt)^y) .{19) 

Further generalizations for including the finite e-effect 
may be deduced as shown above for continuous e-shortest 
path betweenness. 

For one-dimensional convex sets 5* and using 
ai{z,w;x,y;£), Eq. (18) reduces to 



b{x,y) — 2 dfi{z) dfi{w)Q{x - z)Q{y - z) X 

xe{w-x)e{w-y). (20) 

In the limit x -^ y we always have b{x,y]e) -^ 6(x;e). 
We note that b{x, y; e) does not require the condition 
8(e— ||x— y||) as is the case for the corresponding discrete 
estimator (Table lll| . Related generalized concepts of co- 
and group betweenness have been described for discrete 
complex networks |56j . 

For general non-pathological S we almost surely have 
17(2:, w) = 1, i.e., the probability that there are more than 
one global geodesies connecting z and w drawn randomly 
from S according to p is zero. For example, in both 
Figs. [SK and B the set of pairs z, w with a{z, w) = 2 or 
a{z,w) = 4, respectively, is of measure zero. In these 
cases, b{x]e) and b{x,y;e) reduce to 



&(x;e)= // dfi{y)dfi{z)a{y,z;x;e) (21) 



and 



^(a;,2/;£)= // d^j,{z)dfi{w)a{z,w]x,y;e). (22) 



2. Global measures 



Definition 13. The continuous e-average path length 
£(e) = 



d^,{x)d^,{y)^^^ (25) 



e ^ I I dfi{x)dii{y)g{x,y). 



measures the expected geodesic distance in units of e be- 
tween two points X and y drawn randomly according to 



From Eq. (251, the equivalence of this formulation of 



continuous average path length to the intensively stud- 
ied problem in probabilistic geometry [S5] of finding the 
expectation value of the distance between two randomly 
drawn points x,y € S according to the probability dis- 
tribution p becomes evident. Our definitions imply the 
relationship 



C{e) — / diJ,{x)c{x;e) 
's 



(26) 



Definition 14. Similarly, the continuous global e- 
efficiency 



dnix)df^{y){9-^) ) (27) 




d^i{x)d^i{y)g{x,y) ^ 



is the inverse of the expected inverse geodesic distance 
between two points x, y drawn randomly according to p 
measured in units of e. 



Here, we have 



£{£) — ( / dfi{x)e{x; e) 



(28) 



More generally, let As^pirj) be the expectation value of 
a power 77 of the geodesic distance g{x,y) between two 
points x,y € S randomly drawn according to p: 

AsAv) = J J dfi{x) dfi{y)g{x, y)". (29) 

Then continuous e-average path length and global e- 
efhciency may be expressed as 



£(£) = £-lAs,p(l) 



and 



Definition 15. The e-diameter 

V{e) = e^^ sup g{x,y) 

x,y£S 



and the e-radius 



TZ{e) — e ^ ini sup g{x,y) 



(30) 



(31) 



(32) 



(33) 



are global geometric characteristics of the set S that are 
independent of p !l55f . 



E. Further measures 

To illustrate that the proposed framework can be ex- 
tended in several directions, we shortly discuss spectral 
and random-walk-based measures in the context of con- 
tinuous recurrence networks. Motivated by the study of 
eigenvector centrality in complex networks [T], we can 
consider spectral properties of the set S and probability 
density p. 

Definition 16. The linear Laplace operator 



iLJ){x) = / d/i(y)[e(e-||x-y|| 
Js 

-5{x-y)p{y-e)]f{y) (34) 

is a continuous equivalent of the discrete Laplacian ma- 
trix in network theory "[Tj. We are interested in its eigen- 
functions f(x) and eigenvalues A satisfying 



(L,/) {x) = \f{x) 



(35) 



for all X ^ S. 



For example, considering an arbitrary S with uni- 
form p, one obtains an eigenfunction f{x) = C for some 



C G M associated to the eigenvalue A = 0. This is anal- 
ogous to the eigenvector (1,1,...,1) with eigenvalue 
which is always present for the discrete Laplacian matrix 
of general networks [T] . We can expect more interesting 
results for non-uniform p. For example, one may define 
a continuous analog of the eigenvector centrality of com- 
plex network theory 1591 by considering the eigenfunction 
f{x) corresponding to the largest eigenvalue A. 

For discrete networks, there are several measures of 
betweenness based on random walks rather than short- 
est paths [60, 61 . Continuous versions of these measures 
would be based on continuous analogs of random walks 
on 5* that start and end at points y and z randomly cho- 
sen from p. Since in a discrete network the limit distribu- 
tion of a random walk without a sink is proportional to 
the degree distribution, a natural choice for a continuous 
analog is an Ito diffusion process, the limit distribution 
of which is proportional to p{x;e), with a source at y 
and a sink at z J62I . Such a process can most easily be 
defined as a gradient flow dXt = -V^{Xt) dt + \/2T dBt 
that combines a Brownian motion B with a local drift 
coefficient —V'^{Xt) which comes from a potential "^{x) 
that is the product of a temperature T > and the in- 
formation corresponding to p{x]£), which is —hip{x;e). 
The resulting process 



dXt 



T^^dt- 



2TdB, 



(36) 



can then be interpreted as a diffusion that drifts in the 
direction of increasing density. The continuous version of 
Arenas' random walk betweenness (6^ would then be the 
expected density of the process at x when the source and 
sink are drawn from p. Similarly, the continuous version 
of Newman's random walk betweenness |61j would be the 
expected absolute value of the resulting flux density at x 
for a random source and sink. 



F. Behavior under afflne transformations 

All continuous measures defined above are based on 
neighborhood relationships in S and geodesic distances 
between points therein. They are therefore invariant with 
respect to the subclass of affine transformations which 
leaves these properties unchanged, i.e., x — >■ Dx + s 
ioi X G S with D being a combination of rotation and 
isotropic scaling operations and s a translation. This 
is to be understood in the sense that for a measure h, 
h{Dx + s; ae) = h{x; e) holds, where a is the scaling fac- 
tor. The measures considered here are generally not in- 
variant under non-isotropic scaling and shear operations. 



III. DISCRETE ESTIMATORS 

Given the continuous framework defined above, we 
are able to treat the commonly used recurrence network 



TABLE II; A summary of standard unweighted network estimators for tiie continuous geometric properties defined in Sec. Ill] 
For a detailed discussion, see [l] |2l 1141 I63j . SP abbreviates "shortest path" . 



Class 



Name 



Definition 



Comments 



Neighborhood-based 
Local 



Degree 

Degree density 

Local clustering coeff. d = 

Matching index fiij 



Pi ~ N-1 ' 

fei(fci-l) 



C, = iff fci < 2. 



Global 



Edge density 
Transitivity 



n - 1 V^ A-. 



T 



'^i,j,k = l -^ki^kj 



Global clustering coeff. C = jf^^Ci 



Assortativity 



A = 



iE,>.i{*?+fc|)A„-{i(fc,+fc,)) 



^— i = X]i>i ^»j is the number of edges, 

\ 'I.J 



m 



Path-based 
Local 



Global 



Closeness 

Local efficiency 
SP betweenness 






-i j^h^i a 






SP edge betweenness 6,, = C^-^) "' E^,#.,, %lf^ 

Average path length L = (/ij)^ 

Global efficiency E = ( (^Z^) , 

Diameter T) = majcij (Z^j) 

Radius 7?. = min^ maxj (Zij) 



Set lij = A'^ — 1 iff ^ path between i, j j64| 



(Tjft shortest paths connect vertices j, fc, 

o-jk{i) of those include i [M] . 

and Ojk{i,j) include i,j. 

Set Zij = A'^ — 1 iff J path between i,j [T]. 



quantifiers [SI [H] taken from standard complex network 
theory [U [2] as the most straightforward discrete esti- 
mators of the continuous quantities for a finite number 
of observations A*" and finite e. The discrete estimators 
will be denoted using hats, e.g., the discrete estimator 
of continuous average path length C[e) is £(e, N) (we 
will in the following omit the estimators' dependency 
on £ and A^ to simplify the notation). Their numer- 
ical properties have been elaborated in detail in ear- 
lier works US [H [Ml EHl I33I- The characteristics 
of these standard measures for discrete and finite com- 
plex networks have also been studied for random geomet- 
ric graphs and more general network models with strong 
spatial contraints [Sni[S2[53], e.g., the degree distribu- 
tion [50|, network motifs [53], as well as clustering coef- 
ficient and degree correlations [55]. 

Here, we briefly review the estimator's definitions (Ta- 
bleln]). For some specific examples, the estimators will be 
compared to the results theoretically derived from their 
continuous counterparts in Sec. |IV| This will also allow 
us to gain certain insights into their bias and variance for 
finite data sets. 



Weighted network statistics and node splitting 
invariant measures 



We may now ask how the estimation of the above de- 
fined continuous geometric properties from a finite data 
set can be improved with respect to the measures from 
complex network theory that have been used so far for 
this purpose. One way to go in line with standard es- 
timation theory is node-weighted network statistics, as 
proposed by Heitzig et al. [32] • For a full application of 
that theory, weights Wi for all vertices i have to be chosen 
in a suitable way, which we leave as a subject of future 
research. But even with constant weights Wi = 1, the ax- 
iomatic theory developed in [IS] allows us to improve esti- 
mation by using so-called node splitting invariant (n.s.i.) 
versions of network measures to reduce the estimation 
bias that results from excluding self-loops from the net- 
work. Let us illustrate this for the case of continuous 
e-degree density, p{x] e) = /g , x dfi{y). If a; is a vertex, 
p is approximately constant in B^(x), and the latter con- 
tains ki additional vertices (see Table ]ll[), then p{x;e) w 
p{x)Yol{B^{x)) « {ki -\- 1)/N since B^^x) contains k* = 



10 



ki + 1 out of N vertices. In other words, k* /N is a bet- 
ter estimator for p{x; s) than ki/N, since the latter has a 
bias of ~ 0(1 /N). Likewise, the transitivity measure T 
(see Table El) can be improved by using instead the n.s.i. 

transitivity r* = ^^,,.=1 44^^^./ Sf,,fe=i ^^.^^,, 
where A^j — Aij + 6ij = Rij , showing that this approach 
is also more in line with recurrence plot analysis. This 
would also reduce the bias in the estimation of the tran- 
sitivity dimension that was observed in [26, Fig. lOA]. 
The measures k* and T* are examples of n.s.i. measures 
with unit weights, which can basically be interpreted as 
variants of the classical measures in which vertices are 
considered to be linked to themselves HHI. 



B. Choice of the recurrence threshold s 

A careful choice of the recurrence threshold e is crit- 
ical for faithfully estimating the continuous recurrence 
network properties defined above [30]. For too large e, 
i.e., on the order of the diameter of S, boundary effects 
dominate, the discrete recurrence network used for esti- 
mation becomes too dense and is unable to capture the 
geometry induced by S and p (see Sec. |IV A4 ). In con- 
trast, the network's giant component breaks down for 
too small £ with a phase transition at the critical value 
Sc- This obstructs our ability to properly estimate meso- 
scopic and path-based measures for e < Ec- Therefore we 
expect a good performance of the discrete estimators for 
thresholds just above the critical Sc, where much of the 
geometric fine structure is still resolved pS] . 

The problem of selecting e therefore reduces to deriv- 
ing the percolation threshold £c which is directly related 
to the critical edge density pc = p{ec) of the theory of 
random geometric graphs [51] via Eq. m«. pc is linked to 
the commonly studied critical mean degree Zc by 



Pc = 



N -1' 



(37) 



The Erdos-Renyi graph is the simplest random network 
model [P. Since any pair of vertices is linked with 
the same probability p independent of their distance, 
it neglects the effects of spatial embedding. Therefore 
the Erdos-Renyi model is inadequate for describing d- 
dimensional random geometric graphs, and the corre- 
sponding critical mean degree Zc = 1 [I] turns out to 
be too low except for the limiting case d —>■ 00 [35] (see 
Sec. IV A 4 for an example). Taking into account the ef- 



fects of clustering of vertices induced by the spatial em- 
bedding [5S] yields improved analytical bounds on the 
true Zc obtained from numerical simulations [66] . Exact 
analytical results for arbitrary d are not available so far, 
but Dall and Christensen [55] have empirically found the 
scaling law 

Zc{d) ^ Zc{(X)) + Ad-^ (38) 

from extensive numerical simulations, where Zc{oo) = 1, 
7 = 1.74(2) and A = 11.78(5). Inverting pc = p{£c) 



(which is possible as dp{e)/de > in non-pathological 
situations) yields the associated critical threshold 



£c(d) 



Zcjd) 

N - 1 



(39) 



To our best knowledge this is the most useful result avail- 
able so far for our aim of choosing the recurrence thresh- 
old e. However, one should be aware that the results 
of pS] were obtained for the box S — [0, 1]"^ with uni- 
form probability density p which is the most commonly 
studied setting in random geometric graph theory. When 
considering general S and p they may be appropriate as 
a first educated guess for properly selecting e in line with 
the guidelines discussed in iM, 30j. Deriving analytical 
bounds on Zc for such geometries remains an open prob- 
lem. 



IV. EXAMPLES 

We illustrate the above defined continuous geometric 
quantities and their estimators for paradigmatic exam- 
ples by giving closed-form analytical results and relating 
them to numerical evidence from e-recurrence networks 
constructed from time series. The focus will be on exam- 
ples where all quantities of interest can be calculated ei- 
ther analytically or semi-analytically (relying on numer- 
ical evaluation of some integrals), i.e., possessing smooth 
sets S and density functions p{x). This implies that when 
considering the Euclidean norm (which we will use for all 
examples below) and neglecting boundary effects [14m2^ . 
we obtain 



C{x;e) = 1 



dV(d/2) 



2y/^T{{d+l)/2) 



iFi 



1 1-d 3 1 



d + 1 
=C(e) = T{e) 



1 fl-d d + l d + 3 1 

2-ri 



2' 2 ' 2' 4 
+ 3 1 
2 '4 



(40) 



for all transitivity-based measures, where 2^1 (■) is the hy- 
pergeometric function and d the manifold dimension of S. 
A simpler exponential scaling with d can be found for the 
supremum metric. Nontrivial transitivity-based proper- 
ties for fractal sets S and densities p{x) allowing for non- 
integer d, where an analytical calculation of path-based 
measures is problematic, have been treated exhaustively 
in [26.^. 

The results given here hold in the limit £ — ?> 0. 
For simplicity we ignore boundary effects which have 
been treated in [lH [26]. In all examples, we use the 
parametrization ai{y,z\x) to compute continuous e- 
shortest path betweenness. 



11 



A. One-dimensional chaotic maps and stochastic 
processes 

All examples considered in this section are defined on 
convex sets S embedded on the real axis. Therefore, the 
geodesic distance of a;, y e M reduces to g{x, y) = \x — y\. 
Since for one-dimensional S the integral Jg dyp(y) \x — 
y\~^ diverges for all p and all a; G 5", we get f (e) = 
and e{x]e) = oo Vx in all examples of this section. In 
contrast, the corresponding integral always converges for 
non-fractal S with d > 2 and general p. 



1. Bernoulli map / uniformly distributed noise 

The Bernoulli map Xn+i = (2a;„) mod 1 defined on 
the interval S — [0, 1) induces the probability density 
p{x) = 1. This yields 

T{s) = C{e)^C{x;e)^^ 



A , Bernoulli map 

10^ 



Logistic map 



c{x;e) 



1 



2e 



1 - 2a; + 2a;2 
b{x) = 2a;(l-x). 

The same results hold for uniformly distributed noise on 
the interval [0, 1], since S and p{x) are identical to those 
of the Bernoulli map (an exemplary calculation for this 
setting is shown in Appendix ^ . This equality clearly 
illustrates that recurrence network analysis is purely 
geometric and, hence, by design masks out the auto- 
dependency structure of dynamical systems. Stochastic 
and deterministic dynamics can be distinguished when 
embedding techniques are used prior to recurrence net- 
work analysis [T4] . 



2. Gaussian noise 

Considering Gaussian noise with zero mean 
and standard deviation a, i.e., p{x) = 



(l/\/2^ra2) 



TTcr^ exp I 



V(2a2)), 



S = (— oo, +oo) we obtain 

Tie) = Cie)=C{x;e) 
Lie) = ^e- 
c{x;e) = 

b{x) = 



on the real 

3 

" 4 



- a exp I 



1-erf 



"^^'^[vt 




where erf(a;) = -^ L e * dt is the error function. The 
results for mean x ^ can be derived by substituting 




FIG. 4: (Color online) Continuous e-average path length £{e) 
for (A) the Bernoulli map and (B) the logistic map. Analyt- 
ical results are indicated by solid red lines. Estimates jC.{e) 
have been obtained from e-recurrence networks constructed 
from one realization of A = 1,000 samples, respectively, for 
each map (black dots). Ensemble mean (dashed-dotted black 
line) and standard deviation (gray band) for different e have 
been obtained from an ensemble of 100 realizations of each 
model with initial conditions uniformly distributed in the in- 
terval [0, 1]. Vertical dashed lines indicate the estimated per- 
colation thresholds Sc- 



a; —> a: — X on the right side of the equations for the local 
measures given above (see also Sec. II F). 



3. Logistic map 

We can also give exact analytical solutions for the logis- 
tic map in the fully chaotic regime, Xn+i = 4a;„(l — a:„), 
defined on the interval S = [0, 1]. Using the probability 

density p{x) ~ n^^ ^/x{l^^x) [67] yields 



Tie) = Cie)=Cix;e) = - 
Cie) = -^e-' 



zix;e) = 7re[2^a;(l - x) + (1 - 2a;) x 

X (arccos (\/a;) ^ arcsin (v^)) 
8 Im [arcosh (v^)] arcsin (v^) 



b{x) 



4- Comparison to numerical results 



Within an intermediate range of e, the continuous e- 
average path length i2(e) is approximated well by the es- 
timators £(e) calculated from e-recurrence networks for 
both the Bernoulli and logistic map (Fig. [4]) . For small 
e, the estimator breaks down due to the finite number 
of samples used (finite size effect) after the network's gi- 
ant component decomposes into snraller and smaller dis- 
connected components. The Erdos-Renyi approximation 
yields a critical percolation threshold Ec = l/(2(Af— 1)) pa 
5 • 10~^ for both maps using the parameters of Fig. |4J 
which is one order of magnitude smaller than the nu- 
merically observed phase transition point (Fig. l4]). As 



explained in Sec. IIIB this is because the Erdos-Renyi 



12 



model does not account for the effects of spatial em- 
bedding and clustering. To be able to use the relation- 
ship of Eq. (39) based on the empirical results of [2S] 
for 771 = 1, we approximate p(e) — 2e (Bernoulli map) 
and p{e) — 8£artanh(l — 2e)/7r^ (logistic map) for small 
e. This yields Ec ~ 6.4 • 10^"^ for the Bernoulli map 
and £c ~ 6.2 • 10"'^ for the logistic map, which is con- 
sistent with the phase transition points observed numer- 
ically (Fig. |4|. The good agreement of predicted and 
observed phase transition for the Bernoulli map can be 
explained by the fact that the latter exactly meets the as- 
sumptions underlying the theory of Dall and Christensen 
(Sec. 



Bernoulli map 



III B I . These observations indicate that Eq. ( 39 ) is 



indeed useful to derive an educated guess for the proper 
choice of e, even for strongly varying probability densities 
p. Moreover, the phase transition for the logistic map 
occurs at notably larger e than for the Bernoulli map. 
Consistently with the results of ^66], this indicates that 
the increased spatial clustering induced by peaks in the 
density p leads to larger values of the critical mean degree 
Zc and therefore the associated percolation threshold Sc- 
However, the results of [S2] suggest that there may be in 
fact no true phase transition in giant component size for 
non- uniform p in the limit N ^- oo. 



For large £, the approximation in Eq. ( 11 ) is not valid 
anymore and, hence, the discrete estimator breaks down 
in this regime. Note that C{ec') = 1 for a critical Ec'- 
Since £(e) < 1 whereas £{e) = 1 for e > Ec', the defini- 
tion of the discrete estimator is not meaningful anymore 
for thresholds larger than the critical threshold. For the 
Bernoulli map we have Ec/ ~ 1/3 and for the logistic map, 
Ec' = 4/7r^ follows. 

The continuous e-closeness c{x] e) is approximated well 
by the estimator c{x; e) for both the Bernoulli and logistic 
maps (Fig. I5|. However, c(x;e) is notably smaller than 
the true theoretical value particularly in the center of S 
at X — 1/2, implying that shortest paths are longer in 
the empirical e-recurrence network than expected theo- 
retically. This is clearly a finite size effect as the bias and 
variance of the estimator clearly decrease for growing N 
and fixed e (Fig. Isl) . 

The shape of continuous e-shortest path betweenness 
b{x) is approached well by the estimator b{x; e) for both 
maps (Fig. [6]). However, there is a large bias that in- 
creases with the number of samples N, while the variance 
decreases with growing N. That the estimator b{x;E) is 
generally smaller than the theoretical value b{x) for all 
X can be explained by the skipping of vertices due to 
the finite e in the empirical e-recurrence network. This 
effect is expected to increase for growing N when e is 
fixed, since more and more vertices can be skipped along 
a shortest path for the same recurrence radius e, which 
also explains the growing bias in this setting. Accord- 
ingly, the bias decreases for decreasing e when N is suffi- 
ciently large regarding the discussion of suitable choices 
of e in Sec. IHB (Fig. [6]GI,H). However, the bias is not 
a problem in practical situations, because for local mea- 
sures we are usually only interested in relative differences 






Logistic map 
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FIG. 5: (Color online) Continuous e-closeness c{x;e) for 
(A,C,E) the Bernoulli map and (B,D,F) the logistic map. An- 
alytical results are indicated by red solid lines. Estimates 
c{x; e) have been obtained from e-recurrence networks at e = 
0.01 constructed from single reahzations of (A,B) A'^ = 1, 000, 
(C,D) N = 2,000 and (E,F) A^ = 5,000 samples (all: black 
squares). Ensemble mean (black dashed-dotted lines) and 
standard deviation (gray bands) have been calculated as in 
Fig. [4] The standard deviation is too small to be visible in 
the plots for TV = 5, 000 (E,F). 



between vertices and not in the absolute values. 

B. Periodic and quasi-periodic dynamics 

1. Periodic orbit 



We analyze next a periodic orbit (gen- 
eral closed curve) of curve length I embed- 
ded in an m-dimensional phase space, i.e., 
S = {a;eM™:x = /(s);sG[0,/];/(0) = /(/)}, with 
uniform probability density p{x) — \/l. The geodesic 
distance of two points x{s),x{t) along the curve is then 
given by g{x{s), x{t)) = \s — t\. This yields 



r(e) 

£{e) 

c{x;e) 
b{x) 



CiE)^C{x;E)^- 



4:E 

T 
1 

4' 
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FIG. 6: (Color online) Continuous e-shortest path between- 
ness b{x) for (A,C,E,G) the Bernoulli map and (B,D,F,H) the 
logistic map. Analytical results are indicated by red solid 
lines. Estimates b{x; e) have been obtained from e-recurrence 
networks at e = 0.01 constructed from single realizations of 
(A,B) iV = 1,000, (C,D) N = 2,000 and (E,F) N = 5,000 
samples (all: black squares). Panels (G,H) show results for 
£ = 0.005 and A'' — 10, 000. Ensemble mean (black dashed- 
dotted lines) and standard deviation (gray bands) have been 
calculated as in Fig. Ul Note that b{x) is the discrete shortest 
path betweenness of [l4l . but normalized by its theoretical 
maximum value {N - 1)(A - 2)/2. 



As the periodic orbit is a one-dimensional set 
we have £{s) — and e(x;e) = oo as above. 
For example, a circular orbit of radius R as 
generated by a harmonic oscillator with S = 
{x e ]R2 : xi = Ram{s/R),X2 = Rcos{s/R);s £ [0,27ri?]} 
and p{x) ~ 1/{27tR) gives the above results with 
I = 27rR. 



2. Flat 2-torus 

Quasi-periodic dynamics is displayed by a system os- 
cillating with two incommensurable frequencies wi and 
a;2, i.e., where the ratio uji/uj2 is not a rational num- 



ber. The phase space trajectory fills a 2-torus 5* — 
{x = {s,t) : s € [0,2TTR\,t € [0,27rr]} uniformly with 
p{s,t) — p = l/(47r^i?r). The radii R,r are related to 
the oscillation's amplitudes. With the geodesic distance 



',{{s,t),is\t')) 



we obtain 



min (I s — s' 1 , 27ri? — \s — s'\) 
+ min(|t-t'|,27rr-|i-t'|)^)' 



Tie) = C{e)=C{x:e) = l- 



3\/3 

47r 



0.5865 



de) 



ne 



12rR 
-|-2i?^artanh 

— r'^artanh 



ArR\/r^ + R^ + Br''' arsinh 



R 



Vr2 + i?2 

R 



Eie) = 27rrfl£~M 2rartanh ( — 



-i?ln 



r -I- \/r2~+~R2 
-r + Vj-^ + R^ 



Because of symmetry, the local path-based measures do 
not depend on x as for the periodic orbit discussed above 
and we have 

c{x;e) — C{e)^'' 
e{x;e) = ^(e)^"^ 

As expected, in the limit r — > the average path length 
converges to the value obtained for a circle of radius R 
(see above), i.e., liuir^o C{e) — {■KR/2)e~^. 



C. Higher-dimensional symmetric sets 

The ?n-dimensional hyperball and hypercube may be 
viewed as tractable idealizations of higher dimensional 
attracting sets of dynamical systems (in this section, we 
set d = m, since the manifold dimension of S is inte- 
ger). Their study highlights that continuous path-based 
measures may depend sensitively and non-trivially on the 
global geometry of the set. In contrast, their neighbor- 
hood and transitivity-based counterparts just depend on 
the dimension m and are therefore identical for the hy- 



perball and hypercube (Eq. (40)). The sets considered 
here are convex, hence, g{x,y) — \\x — y\\2 holds when 
using the Euclidean norm. 



1. m-dimensional hyperball 

Here we consider the m-dimensional hyperball S 
'" with the uniform probability density p{x) = p 
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£L(6) 




FIG. 7: Continuous e-average path length eC{e) of the hy- 
perbaU S'" with uniform probability density, obtained from 



Eq. (421 



1/ Sgdx = Vol(S™) ^ == r(f + l)/7rt. Following 
Hanimersley [68], the r-th moment of the distribution 
of point-to-point distances ||a; — 2/II2 in §'" is given by 



^^n 



r (im + i) (m + r)T {m + \r + l) 
Then the continuous e-average path length is 



(41) 



£(e) 



lie 



r(TO + i) r(iTO + i) 
lyr(m-fl)' 



ir(im 



(42) 



and all its higher moments are known via Eq. ( |41[ ) (see 
Fig. I7|. Some examples for lower-dimensional spheres 
follow: 



eC{e) 



-IT 



= < 



2 

3 

128 , 

45 
36 
35 
16384 , 

4725 ' 
800 
693 



771 = 1 
m — 2 

m = S 
4 
5 



m 



m 



Note that the result for to = 1 agrees with the cor- 
responding one for the Bernoulli map when consider- 
ing the stretching of the domain by a factor of 2, since 
S"'^ = [—1,1]. In the limit to, — > 00 the continuous e- 
average path length is C{e) — -\/2e^^ (see Fig. [7J). We 
can also derive in closed form an expression for the con- 
tinuous £-closeness c(0; e) of the center of §™, taking ad- 
vantage of the spherical symmetry: 



c(0; e) ^ — e ^ dxi . . . dxm P\ x\ 



= e-^VL„,p / drr'^'^r. 



With the full solid angle in to dimensions i7„ 

77i7rT/r(f + 1) this leads to 

c(0;e) — e. 



(43) 



Note that the limit £c(0;e)~^ — >■ 1 for tti — ?► 00 shows 
that almost all of the measure ^(S™) of the unit radius 
hyperball S™ is concentrated at its surface for large m. 
For the special case of to = 2 (unit disk with uniform 
p(x)), Lew et al. |69j give a nearly closed-form expression 
for the continuous e-closeness at x{q) — {q,0), 

c{x{q);er' - ^ (I6(g2 - l)K{q') + 4{q' + 7)E{q^)) £-\ 

where < g < 1 and the value for arbitrary a; e S^ 
may be obtained after an appropriate rotation. K{ni) 
and E{m) are complete elliptic integrals of the first and 
second kind (see §17.3 in 70J). 

For the continuous local e-efficiency e(0; e) of the center 
of S™ we get for to > 1 



e(0;e) — £ dxi . . .dxmPy x\ 



= eVlraP / drr'^~^r~^ 
Jo 

TO, 
771—1 

A somewhat more involved calculation of the continuous 
e-betweenness 6(0) of the center of S™ yields (see Ap- 
pendix Icj) 



6(0) = 



1 



r(f + 1) 

777 TTT 



(44) 



The high degree of symmetry of S™ allows us to derive 
closed-form results for local path-based measures at its 
center for many p{x) = p{r, fl), as long as the probability 
density separates into a radial and an angular part, i.e., 
p{r,n) =p{r)p{n). 



2. m-dimensional hypercube 



The hypercube S = K™ = [0, 1]™ with uniform prob- 
ability density p{x) = Vol(K™)^^ = 1 is much harder to 
treat analytically than the hyperball S™. Hence, rigorous 
results are only available for isolated dimensions 777 and a 
subset of the continuous measures defined above [5H1 FfTI - 
I73j . Solving the resulting general box integrals remains 
a largely unsolved problem of applied and experimental 
mathematics. 

The following closed-form expressions for the continu- 
ous £-average path length £(e) are based on the expecta- 
tion values for point-to-point distances Ak>7i(1) — e£{e) 
(see Eq. (l29|) listed in [7S]: 



eC{e) 
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TO : 


= 1 


ji; (2 + 72 + 5 In (1 + 72)) 


TO : 
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_118 _ 2 , 34 ./2 - 
21 3"^ ^ 21 V^ 


-|V3 


TO = 
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+21n(l + \/2)+8 


M^) 
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FIG. 8; (Color online) Continuous e-average path length 
eC{e) of the hypercube K™ with uniform probability density, 
obtained by numerical Monte Carlo integration using Mathe- 
matica (yellow diamonds). Analytical lower (blue discs) and 
upper (red squares) bounds from ^T^ are also shown. 



Some numerical results for ni = 1, . . . , 10 are displayed 
in Fig. Is] Andcrssen et al. 71j proved the bounds 



<eC< 



'\ 



1+2W1- 



5m 



(45) 



implying e>C — > oo for to — > cxd. This is in contrast to the 
hyperball, where this limit is finite (see above) . Using ex- 
pectation values for the inverse point-to-point distances 
Ak™(— 1) — e^^£{e)^^, we are able to give the following 
expressions for the continuous e-efficiency £{s): 



{e£{s))-' = { 



1(1 - V2) + 4\n{l + ^ 


i-i- + i^/2-|V3 


+21n(l + \/2) 


+121n(i±#) 


-41n(2 + \/3) 



V2) 



TO, = 1 



m = 3 



Note that as 5 = K-'^ = [0,1], the results for to, = 1 
agree with the corresponding ones for the Bernoulli map 
for both continuous e-average path length and efficiency. 
Further expressions for Ak.t.(1) and AK'n(— 1) for to = 
4, 5 are given in [73]. 

Another object of interest in the theory of box integrals 
is the integral 



S™(7/) 



dxllxii" 



2' 



(46) 



which is related to the continuous e-closeness c(0; e) of 
the origin a; = (and, by symmetry, to that of all the 
2™ corners of the hypercube) for 77 = 1 and to the local 
efficiency e(0;£) of the same points for 77 = —1: 



c{0;e) = B„-,il)-'e 
e(0;e) = B^{-l)e. 



(47) 
(48) 



We can now once again use results from Bailey et al. 
to give some closed forms for small m: 



(rc(0;e)-^ 



( 1 

2 



hV2 + Hi + V2)) 



3 

IV3 



m = 1 
TO = 2 



24' 



i ln(2 + ^/3) m = 3 



and 



-'e{0;e) 



00 

21n(l + %/2) 



m 
m 



1 

2 . 



-\- 



ln(2 + V3) TO = 3 



Further solutions for to = 4, 5 are given in [73] . 



V. DISCUSSION 

We have shown that the definitions of continuous 
geometric measures provided in this paper are feasi- 
ble for describing e-recurrence networks for time series 
analysis as well as, more generally, random geometric 
graphs [551 [SOI [51] ■ Our theoretical framework may 
readily be generalized to encompass other classes of ran- 
dom networks with spatial constraints (i.e., spatial net- 
works) [ini [511 US] J e.g., those with an edge length distri- 
bution of the form pi{l) oc exp(— Z/^) describing among 
others the substrate of climate networks [22 [^ (in con- 
trast to the sharp cutoff pi(Z) ex Q{e — I) for random ge- 
ometric graphs). For the neighborhood-based measures, 
this generalization can be achieved by substituting terms 
containing the Heaviside function or B^ [x) with suitably 
chosen expressions involving pi{l). One possible appli- 
cation to real-world spatial networks is computing ex- 
pectation values for the characteristics of an ensemble of 
spatial random network surrogates to assess which prop- 
erties of a given empirical network can be explained by 
Pi{l) alone. Additionally, more general metrics could be 
used for measuring the distance I between connected ver- 
tices. Research along these lines may also help to shed 
light on the specific topology and dynamics of growing 
spatial complex networks (cf. [HI [75] ) . 

Furthermore, we have demonstrated that the result- 
ing continuous properties can be approximated by esti- 
mators calculated from empirical £-recurrence networks 
reasonably well, even for relatively small N and large 
e. The continuous framework promotes considerable ad- 
vances in the theoretical understanding of e-recurrence- 
network-based time series analysis. Among others, from 
the examples of hyperballs and hypercubes in various di- 
mensions TO, the claim that path-based measures depend 
explicitly on the global geometry of the set 5* is theo- 
retically justified. This is in contrast to the continuous 
notions of local and global transitivity, as at least the 
continuous local e-clustering coefficient C{x;s) depends 
on the local dimensionality of the set S [26]. Along 
these lines, in the future we may gain an understand- 
ing of the differing performance of transitivity-based and 
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path-based measures in classifying qualitatively different 
behavior of dynamical systems |131 I33j . For example, 
more complex dynamical systems such as the Lorenz and 
Rossler models or noisy dynamical systems, where no 
closed-form expression for the invariant probability den- 
sity p exists, may be studied by estimating p from simu- 
lated trajectories, e-recurrence network measures could 
then be calculated by numerical integration techniques 
relying on p and the integral expressions given in this 
paper. Circumventing the computational limitations of 
discrete e-recurrence network analysis when N ^ oo, this 
approach would in principle allow us to approximate the 
continuous geometric quantities defined above as closely 
as desired. 

The examples of hyperballs and hypercubes establish 
links to some current research problems in probabilis- 
tic geometry and applied mathematics, among others, to 
the theory of box integrals [58l [711 - F73] . Perhaps these 
highly symmetric model sets could serve to understand 
theoretically some qualitative features of path-based e- 
recurrence network measures for strange attractors such 
as the Lorenz or Rossler attractors. It remains an open 
question as to whether it is possible to solve the inte- 
grals for continuous path-based measures in the case of 
self-similar sets S and more complex, potentially also self- 
similar densities p. 

The theoretical framework put forward in this paper 
enables several practical advances, which are particu- 
larly relevant for applications to time series analysis of 
real-world data. Analytical solutions for continuous e- 
recurrence network measures allow us to assess the bias 
and variance of the discrete estimators from complex net- 
work theory that have been used in the literature so 
far. These insights led to devising improved discrete es- 
timators based on the concept of node-weighted network 
statistics [IS]. Furthermore, for the first time we were 
able to formulate a theoretically motivated criterion for 
the selection of the recurrence threshold e based on the 
critical percolation threshold Ec which for a given system 
can be estimated using our theory. 

Finally, we should note that we now have a compre- 
hensive continuous theory for essentially all relevant mea- 
sures of £-recurrence networks. This foundation will 
help to further increase our understanding as well as 
strengthen the general confidence in the method of e- 
recurrence network analysis in practical situations, e.g., 
the analysis of real-world time series. Our results sug- 
gest that e-recurrence network analysis is the simplest 
and best understood network-based approach to nonlin- 
ear time series analysis available so far. 



,Wk € S with 
u'i — 2:^ 1 1 < S for i — I . . . k — 1, 
Wi\\ < e for i — 1 . . . fc, implying 



cause S is path-connected, l{x, y; e) is finite. We also put 
l{x,x;e) = 0. Let lij{e,N) > l{xi,Xj;e) be the network 
distance between Xi and Xj in the e-recurrence network 
constructed from the first N points of a sequence of in- 
dependent draws from p. One can then prove that for 
fixed £, fixed nodes Xi,Xj, and A^ -> oo, it has probabil- 
ity one that lij{e,N) = l(xi,Xj;e) eventually (i.e., there 
is some N{i,j,e) so that lij{e,N) = l{xi,Xj]e) for all 
iV > N{i,j,e)). This is because for k = l{xi,Xj;e), there 
is 5 > and zq, . . . , Zk € S with zq = x, Zk = y, and 
I |zi_i — Zil I < e ~ 26 for i = 1 . . . k, and with probability 
one, the sequence contains points wq 
Wf) — Xi, Wk = Xj, and 
so that also ||wi_i 
lij{e,N) < k when N > N{i,j,e) where N{i,j,£) is the 
index of the last of the Wi to occur in the sequence. 

Moreover, l{x,y;e) < g{x,y)/e + 1 and el{x,y;e) < 
g{x, y) if g{x, y) is not an integer multiple of e. This is 
because for all J > and e' < e, there is a path from x to 
y of length < g{x, y) + 5, hence for k — \{g{x, y) + 5)/e'^ 
{\x'\ is the smallest integer not less than x), there are 
zq, . . . ,Zfc e S with zq ^ X, Zk = y, and ||zj_i - Zi\\ < 
e' < £ for i — 1 . . . fc, so that l{x,y;e) < fc. On the 
other hand, if S is sufficiently well-behaved, one will also 
have el{x,y;e) /* g{x,y) for e — > 0. More precisely, 
assume S is "locally almost convex" in the sense that for 
all i > 1, there is some e > so that for all x,y S S 
with ||a; — j/ll < e, we have g{x,y) < Le. Then for all 
L > 1, there is some e > so that el{x, y; e) > g{x, y)/L. 
Putting all these facts together, we see that elij{e,N) is 
a plausible estimate oi g{x,y). 



Appendix B: Continuous e-average path length for 
Bernoulli map and uniformly distributed noise 

For illustration, we give the detailed calculation of£{e) 
for the Bernoulli map and, equivalently, uniformly dis- 
tributed noise: 



r£(e) = / / dxdy\x~y\ 
Jo Jo 



dx 



dx 



dx 



dy\x-y\+ / dy \x - y\ 







dy{y-x)+ / dy {x - y) 
Jo 
-i 1 



,y "Xy 



1 2 

xy- ^y 



0, 



Appendix A: Sketch of proof of Eq. ( 11 ) 



For e > 0, we define the continuous s-distance l(x, y; e) 
between x ^ y G S to he the smallest integer fc > 
such that there are points zg, . . . , zj, € S with zq = x, 
Zk = y, and \\zi-i — Zi\\ < e for z = 1 . . . fc. Note that he- 



dx 



dx 



1 1 2 , 2 I 2 1 2 

- - .T - -X +X +X - -X 



X + X^ 



1 3 1 2_^1 

-X — -rX H- -X 



n 1 
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Appendix C: Continuous e-betweenness for tlie 
center of a hyperball 
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